Partial LD 1
Hongtong
2024-11-09
R Markdown
This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see http://rmarkdown.rstudio.com.
When you click the Knit button a document will be generated that includes both content as well as the output of any embedded R code chunks within the document. You can embed an R code chunk like this:
library(ggplot2)
library(ggprism)
library(tidyverse)
library(dplyr)
library(plotly)
library(EnvStats)
library(car)Plotting examples
#read the LD csv file
(partial_LD <- read_csv("/Users/cynthialin/Desktop/School/lab/Chao Lab/Atoh1-Ebf3 behavior/20240808_Partial LD.csv"))## # A tibble: 144 × 25
## Cohort Batch Subject `Home cage` ID Sex Genotype Zone
## <chr> <chr> <chr> <chr> <chr> <chr> <chr> <chr>
## 1 C6 Batch 4 Cage 1 BEAF18 BEAF18-1 F C Light
## 2 C6 Batch 4 Cage 1 BEAF18 BEAF18-1 F C Dark
## 3 C6 Batch 4 Cage 2 BEAF18 BEAF18-2 F D Dark
## 4 C6 Batch 4 Cage 2 BEAF18 BEAF18-2 F D Light
## 5 C6 Batch 4 Cage 3 BEAF18 BEAF18-3 F E Light
## 6 C6 Batch 4 Cage 3 BEAF18 BEAF18-3 F E Dark
## 7 C6 Batch 4 Cage 4 BEAF18 BEAF18-4 F B Light
## 8 C6 Batch 4 Cage 4 BEAF18 BEAF18-4 F B Dark
## 9 C6 Batch 5 Cage 1 BEAF19 BEAF19-1 F F Light
## 10 C6 Batch 5 Cage 1 BEAF19 BEAF19-1 F F Dark
## # ℹ 134 more rows
## # ℹ 17 more variables: `Duration (Crossover)` <dbl>,
## # `Entry Count (Crossover)` <dbl>, `Latency (Crossover)` <dbl>,
## # `Total Distance (Crossover)` <dbl>,
## # `Horizontal Activity Count (Crossover)` <dbl>,
## # `Ambulatory Activity Count (Crossover)` <dbl>,
## # `Rest Time (Crossover)` <dbl>, `Rest Episode Count (Crossover)` <dbl>, …show the first part of the csv file
head(partial_LD)## # A tibble: 6 × 25
## Cohort Batch Subject `Home cage` ID Sex Genotype Zone
## <chr> <chr> <chr> <chr> <chr> <chr> <chr> <chr>
## 1 C6 Batch 4 Cage 1 BEAF18 BEAF18-1 F C Light
## 2 C6 Batch 4 Cage 1 BEAF18 BEAF18-1 F C Dark
## 3 C6 Batch 4 Cage 2 BEAF18 BEAF18-2 F D Dark
## 4 C6 Batch 4 Cage 2 BEAF18 BEAF18-2 F D Light
## 5 C6 Batch 4 Cage 3 BEAF18 BEAF18-3 F E Light
## 6 C6 Batch 4 Cage 3 BEAF18 BEAF18-3 F E Dark
## # ℹ 17 more variables: `Duration (Crossover)` <dbl>,
## # `Entry Count (Crossover)` <dbl>, `Latency (Crossover)` <dbl>,
## # `Total Distance (Crossover)` <dbl>,
## # `Horizontal Activity Count (Crossover)` <dbl>,
## # `Ambulatory Activity Count (Crossover)` <dbl>,
## # `Rest Time (Crossover)` <dbl>, `Rest Episode Count (Crossover)` <dbl>,
## # `Movement Time (Crossover)` <dbl>, …show the variable names of the csv file
names(partial_LD)## [1] "Cohort"
## [2] "Batch"
## [3] "Subject"
## [4] "Home cage"
## [5] "ID"
## [6] "Sex"
## [7] "Genotype"
## [8] "Zone"
## [9] "Duration (Crossover)"
## [10] "Entry Count (Crossover)"
## [11] "Latency (Crossover)"
## [12] "Total Distance (Crossover)"
## [13] "Horizontal Activity Count (Crossover)"
## [14] "Ambulatory Activity Count (Crossover)"
## [15] "Rest Time (Crossover)"
## [16] "Rest Episode Count (Crossover)"
## [17] "Movement Time (Crossover)"
## [18] "Movement Episode Count (Crossover)"
## [19] "Stereotypy Time (Crossover)"
## [20] "Stereotypy Episode Count (Crossover)"
## [21] "Stereotypy Activity Count (Crossover)"
## [22] "Vertical Activity Count (Crossover)"
## [23] "Vertical Activity Time (Crossover)"
## [24] "Clockwise Revolutions (Crossover)"
## [25] "Counter-Clockwise Revolutions (Crossover)"(lightPhase <-
partial_LD %>%
filter(partial_LD$Zone == "Light"))## # A tibble: 72 × 25
## Cohort Batch Subject `Home cage` ID Sex Genotype Zone
## <chr> <chr> <chr> <chr> <chr> <chr> <chr> <chr>
## 1 C6 Batch 4 Cage 1 BEAF18 BEAF18-1 F C Light
## 2 C6 Batch 4 Cage 2 BEAF18 BEAF18-2 F D Light
## 3 C6 Batch 4 Cage 3 BEAF18 BEAF18-3 F E Light
## 4 C6 Batch 4 Cage 4 BEAF18 BEAF18-4 F B Light
## 5 C6 Batch 5 Cage 1 BEAF19 BEAF19-1 F F Light
## 6 C6 Batch 5 Cage 2 BEAF19 BEAF19-2 F D Light
## 7 C6 Batch 5 Cage 3 BEAF19 BEAF19-3 F C Light
## 8 C6 Batch 5 Cage 4 BEAF19 BEAF19-4 F E Light
## 9 C6 Batch 6 Cage 1 BEAF20 BEAF20-1 F B Light
## 10 C6 Batch 6 Cage 2 BEAF20 BEAF20-2 F F Light
## # ℹ 62 more rows
## # ℹ 17 more variables: `Duration (Crossover)` <dbl>,
## # `Entry Count (Crossover)` <dbl>, `Latency (Crossover)` <dbl>,
## # `Total Distance (Crossover)` <dbl>,
## # `Horizontal Activity Count (Crossover)` <dbl>,
## # `Ambulatory Activity Count (Crossover)` <dbl>,
## # `Rest Time (Crossover)` <dbl>, `Rest Episode Count (Crossover)` <dbl>, …mod1 <- lm(`Duration (Crossover)`~Sex+Genotype+Sex*Genotype, data = lightPhase)
(car::Anova(mod1,type = 2 ))## Anova Table (Type II tests)
##
## Response: Duration (Crossover)
## Sum Sq Df F value Pr(>F)
## Sex 1404 1 0.2290 0.6340
## Genotype 22828 5 0.7447 0.5931
## Sex:Genotype 9194 5 0.2999 0.9110
## Residuals 367840 60par(mfrow = c(2, 2))
plot(mod1)
anv <- aov(`Duration (Crossover)`~Sex+Genotype, data = lightPhase)
(posthoc <- TukeyHSD(x = anv, conf.level=0.95))## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = `Duration (Crossover)` ~ Sex + Genotype, data = lightPhase)
##
## $Sex
## diff lwr upr p adj
## M-F -8.317778 -44.16913 27.53358 0.6446608
##
## $Genotype
## diff lwr upr p adj
## B-A 30.941237 -62.41227 124.29475 0.9247653
## C-A 15.575404 -77.77811 108.92891 0.9963715
## D-A 12.381515 -79.23870 104.00173 0.9986721
## E-A 2.856729 -88.76349 94.47694 0.9999991
## F-A -28.647980 -124.00932 66.71336 0.9494770
## C-B -15.365833 -106.66737 75.93570 0.9962198
## D-B -18.559722 -108.08824 70.96880 0.9900395
## E-B -28.084509 -117.61303 61.44401 0.9396768
## F-B -59.589217 -152.94273 33.76429 0.4269340
## D-C -3.193889 -92.72241 86.33463 0.9999981
## E-C -12.718675 -102.24720 76.80985 0.9983119
## F-C -44.223384 -137.57689 49.13013 0.7321755
## E-D -9.524786 -97.24446 78.19489 0.9995404
## F-D -41.029495 -132.64971 50.59072 0.7757400
## F-E -31.504709 -123.12492 60.11551 0.9130692#base ggplot
plot <- lightPhase %>%
mutate(Genotype = factor(Genotype)) %>%
ggplot2::ggplot(
aes(x = Genotype, y = `Duration (Crossover)`, fill = Genotype)
)+
list(
geom_violin(width = 0.5),
geom_boxplot(width = 0.1, color = "black", alpha = 0.2),
geom_jitter(aes(shape = lightPhase$Sex), width = 0.2), # jitter categorized by sex
scale_shape_manual(values = c(1, 16)), #1 is open circle; 16 is closed circle
stat_n_text(), # sample size calculation
theme(
legend.position = "right",
plot.title = element_text(size=11)
),
ggtitle("plot"),
xlab("")
)
#base ggplot with prism appearance
plot +
scale_color_prism("floral") +
scale_fill_prism("floral") +
guides(y = "prism_offset_minor") +
theme_prism(base_size = 16) +
theme(legend.position = "right")## Warning: The S3 guide system was deprecated in ggplot2 3.5.0.
## ℹ It has been replaced by a ggproto system that can be extended.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
#interactive plot
idInfo <- paste(lightPhase$ID)
genotypeInfo <- paste("Genotype: ", lightPhase$Genotype)
sexInfo <- paste("Sex: ", lightPhase$Sex)
durationInfo <- paste("Crosscover Duration: ", lightPhase$`Duration (Crossover)`)
plot2 <- lightPhase %>%
mutate(Genotype = factor(Genotype)) %>%
ggplot2::ggplot(
aes(x = Genotype, y = `Duration (Crossover)`, fill = Genotype,
text = paste(idInfo, genotypeInfo, sexInfo, durationInfo, sep = "\n")
)
)+
list(
#geom_violin(width = 0.5),
geom_boxplot(width = 0.1, color = "black", alpha = 0.2),
geom_jitter(aes(shape = lightPhase$Sex), width = 0.2), # jitter categorized by sex
scale_shape_manual(values = c(1, 16)), #1 is open circle; 16 is closed circle
stat_n_text(), # sample size calculation
theme(
legend.position = "right",
plot.title = element_text(size=11)
),
ggtitle("plot"),
xlab("")
)
ggplotly(plot2, tooltip = "text")